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A  reduced  ID  dynamic  model  of  a  planar  direct  internal  reforming  SOFC  (DIR-SOFC)  is  presented  in  this 
paper  for  system  research  by  introducing  two  simplifications.  The  two  simplification  strategies  employed 
are  called  Integration  and  Average,  respectively.  The  present  model  is  evaluated  with  a  detailed  ID  SOFC 
model,  which  does  not  introduce  the  two  simplifications,  and  a  lumped  parameter  (i.e.,  OD)  SOFC  model. 
Results  show  that  under  the  operating  conditions  investigated  the  accuracy  of  the  reduced  model  is 
not  significantly  compromised  by  the  two  simplifications  in  prediction  of  the  outlet  gas  flow  rates  and 
molar  fractions,  the  outlet  temperatures,  and  the  cell  voltage,  while  its  computational  time  is  significantly 
decreased  by  them.  Moreover,  it  is  quite  simple  in  form.  Therefore,  the  reduced  SOFC  model  is  attractive 
for  system  research.  Compared  with  the  lumped  model,  the  reduced  SOFC  model  is  an  improvement 
with  regard  to  accuracy  because  it  takes  into  account  the  spatially  distributed  nature  of  SOFCs  to  a  cer¬ 
tain  extent.  The  discretized  node  number  for  solving  the  reduced  model  can  be  taken  as  an  adjustable 
parameter  in  modeling,  and  is  determined  according  to  specific  modeling  requirements. 

©  2008  Elsevier  B.V.  All  rights  reserved. 


1.  Introduction 

The  solid  oxide  fuel  cell  (SOFC)  is  attracting  more  and  more  inter¬ 
ests  as  a  clean  and  highly  efficient  device  for  power  generation  [1  ]. 
Since  it  uses  an  oxide  ceramic  material  as  the  electrolyte,  high  oper¬ 
ating  temperature  (typically  600-1000  °C)  is  required  to  achieve 
sufficient  ionic  conductivity.  The  high  temperature  makes  it  possi¬ 
ble  that  the  fuel  reforming,  when  a  hydrocarbon  is  used  as  the  fuel, 
is  carried  out  directly  on  the  anode,  giving  rise  to  a  direct  internal 
reforming  SOFC  (DIR-SOFC).  Moreover,  the  high  operating  temper¬ 
ature  also  produces  high  quality  waste  heat,  which  can  be  recycled 
by  heat  recovery  devices  or  gas  turbines  to  make  combined  heat  and 
power  (CHP)  or  hybrid  systems  with  very  high  overall  efficiency. 

In  SOFC  system  research,  simple  and  accurate  SOFC  dynamic 
models  are  very  valuable,  because  the  model-based  simulation  is 
a  valid  tool  for  optimizing  the  system  parameters  and  flowsheet, 
as  well  as  for  developing  control  strategies.  In  recent  years,  a  num¬ 
ber  of  distributed  parameter  SOFC  models,  ranging  from  ID  to  3D 
models,  have  been  presented  [2-4].  These  models  consider  the 
spatially  distributed  nature  of  SOFC  variables  and  the  complex 
transport  processes  within  the  SOFC,  and  are  useful  for  optimiz¬ 
ing  cell  design  and  operating  conditions.  However,  for  SOFC  system 
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research  such  models  have  two  disadvantages.  First,  these  mod¬ 
els  are  quite  complex  in  form,  and  consequently  are  difficult  to 
build  and  use.  Second,  solving  these  models  using  the  numeri¬ 
cal  method  is  usually  time-consuming.  Up  to  now,  what  widely 
used  in  SOFC  system  research  is  the  lumped  parameter  model  (i.e., 
the  OD  model)  [5-7],  which  is  simple  in  form  and  has  an  obvi¬ 
ous  advantage  in  computational  efficiency  because  the  internal 
transport  processes  and  the  spatially  distributed  nature  of  SOFCs 
are  neglected.  However,  the  accuracy  of  the  lumped  parameter 
model  can  be  impaired  for  the  same  reasons.  This  is  particularly 
true  for  the  DIR-SOFC,  which  incorporates  the  reforming  process 
into  the  cell  channels.  Correspondingly,  the  accuracy  and  cred¬ 
itability  of  studies  by  use  of  the  lumped  parameter  model  are  also 
affected. 

This  paper  will  present  a  reduced  1 D  dynamic  model  of  a  planar 
DIR-SOFC  for  system  research,  which  is  simpler  than  the  general  1 D 
model  because  of  introducing  two  simplifications.  Compared  with 
the  lumped  parameter  model,  the  present  model  considers  the  spa¬ 
tially  distributed  nature  of  SOFCs  to  a  certain  extent,  and  thus  is  an 
improvement  with  regard  to  accuracy.  A  detailed  ID  SOFC  model, 
which  does  not  introduce  the  two  simplifications,  and  a  lumped 
parameter  SOFC  model  are  employed  to  evaluate  the  present  model. 
The  rest  of  this  paper  is  organized  as  follows.  In  Section  2,  the 
reduced  SOFC  model  is  presented,  followed  by  the  solution  method 
given  in  Section  3.  In  Section  4,  the  model  evaluation  is  presented. 
Finally,  conclusions  are  drawn  in  Section  5. 


Y.-W.  Kang  et  al.  /  Journal  of  Power  Sources  188  (2009)  170-176 


171 


Nomenclature 

at  k  stoichiometric  vector  of  species  i  in  stream  k  in  reac¬ 
tions  (i)-(iv) 

ak  total  stoichiometric  vector  of  stream  k  in  reactions 

(i)-(iv) 

A  area(m2) 

cvi  molar  specific  heat  capacity  of  gas  species  i 
(J  mol-1  K_1 ) 

cPt pen.  cp  int  mass  specific  heat  capacity  of  the  PEN  and  inter¬ 
connect  (Jkg-1  K-1 ) 

C  gas  molar  concentration  (mol  m-3 ) 

Ea  activation  energy  (J  mol-1 ) 

F  Faraday  constant  (=96,485  C  mol-1 ) 

h  gaseous  molar  specific  enthalpy  (J  mol-1 ) 

AH  reaction  enthalpy  change  at  the  temperature  of 

298.15  K(J  mol-1) 

J  average  current  density  (A  m-2 ) 

kr  steam  reforming  reaction  constant 

(=4274  mol  s_1  m-2  bar-1 ) 

/<sh  WGS  reaction  constant  (=1000) 

I<sh  equilibrium  constant  of  WGS  reaction 

I  cell  length  (m) 

n  node  number 

ne  number  of  electrons  transferred  per  electrochemical 

reaction 

N  molar  flow  rate  (mol  s_1 ) 

P  pressure  (bar) 

Qcond  thermal  power  due  to  heat  conduction  (W) 
rj  reaction  rate  of  reaction  j  (mol  s_1  nrr2) 

R  universal  gas  constant  (=8.314J  mol-1  K-1 ) 

P0hm  cell  internal  resistance  ( £2  m2 ) 

Sm  molar  source  term  (mol  s_1  m-3 ) 

Sh  heat  source  term  ( W  nrr3 ) 

T  temperature  (I<) 

U  cell  operating  voltage  (V) 

H0cv  open  circuit  voltage  (V) 

U0  open  circuit  voltage  at  the  standard  pressure  (V) 

W  cell  width  (m) 

xik  molar  fraction  of  gas  species  i  in  stream  k 
z  spatial  coordinate  (m) 

Greek  symbols 

tjk  number  of  gas  species  in  stream  k 

X  thermal  conductivity  (W  m_1  K-1 ) 

p  density  (kg  m-3) 

Ar  computational  time  per  time  step  (s) 


2.  Model  development 

Fig.  1  shows  the  schematic  diagram  of  a  unit  planar  DIR-SOFC. 
The  unit  cell  consists  of  two  interconnect  plates  and  a  tri-layer 
structure  composed  of  two  porous  electrodes,  anode  and  cathode, 
separated  by  a  dense  ion-conducting  electrolyte  (often  referred  to 
as  the  PEN).  To  compose  an  SOFC  stack,  repeating  unit  cells  need 
to  be  further  connected  in  series.  The  methane  mixed  with  water 
steam  (or  their  product  after  pre-reforming  of  some  degree)  is  usu¬ 
ally  used  as  fuel  at  the  anode  side,  and  air  as  oxidant  at  the  cathode 
side.  Due  to  the  high  temperature  in  the  SOFC,  the  steam  reform¬ 
ing  process  directly  occurs  on  the  anode.  Through  electrochemical 
reactions  the  electricity  is  generated.  All  the  reactions  taken  into 
account  in  this  paper  are  listed  in  Table  1. 

For  the  model  development,  several  general  assumptions  are 
made  as  follows: 

(1)  Each  unit  cell  in  the  SOFC  stack  operates  identically,  and  there 
is  no  heat  transfer  between  adjacent  unit  cells.  Therefore,  a  unit 
cell  can  be  investigated  to  represent  the  full  stack  performance. 

(2)  The  SOFC  is  treated  as  a  ID  plant  along  the  gas  flow  direction. 

(3)  All  gases  are  assumed  to  be  ideal  gas. 

(4)  The  pressure  in  the  gas  channels  is  assumed  to  be  constant. 

To  develop  the  reduced  SOFC  model,  two  simplifications  are 
further  introduced: 

(a)  The  PEN,  interconnects  and  gas  channels  are  integrated  together 
along  the  perpendicular  direction,  i.e.,  the  SOFC  is  considered 
to  have  only  one  temperature  layer. 

(b)  The  current  density  distribution  is  considered  to  be  uniform 
within  the  SOFC,  and  the  cell  voltage  is  determined  by  the  aver¬ 
age  gas  molar  fractions  and  cell  temperature. 

These  two  simplification  strategies  can  be  called  Integration  and 
Average,  respectively.  By  introducing  Simplification  (a),  the  com¬ 
plex  heat  transfer  and  gas  diffusion  between  the  SOLID  (including 
both  the  PEN  and  interconnects)  and  the  gas  phase  are  neglected, 
and  thereby  the  SOFC  model  can  be  greatly  simplified  in  form.  As 
both  the  electrodes  and  interconnects  are  normally  good  conduc¬ 
tors,  the  current  density  actually  is  spatially  distributed  to  ensure 
a  constant  operating  voltage  throughout  the  cell.  To  calculate  the 
current  density  distribution  and  cell  voltage,  the  time-consuming 
iterative  computation  is  needed.  By  introducing  Simplification  (b), 
the  cell  voltage  can  be  directly  calculated  by  the  average  current 
density,  and  the  iterative  computation  is  thereby  avoided.  Treat¬ 
ments  similar  to  Simplification  (b)  can  also  be  found  in  Refs.  [8,9]. 


Subscripts 
a  anode 

c  cathode 

int  interconnect 

i  gas  species  i 

j  reaction  j 

k  gas  stream  k 

PEN  PEN  structure 

r  steam  reforming  reaction 

s  SOLID  (including  both  the  PEN  and  interconnects) 

sh  WGS  reaction 

Superscripts 
in  inlet 

out  outlet 

ref  reference  condition  (298.15  K,  1  bar) 

spatially  average  value 


O 


Fuel 


|  Interconnect 

CH4  +H20->CO  +  3H2  (steam  reforming) 
CO  +  H20  C02  +  H2  (water-gas  shift) 


Depleted 

fuel 


Fig.  1.  Schematic  diagram  of  a  unit  planar  DIR-SOFC. 


172 


Y.-W.  Kang  et  al.  /  Journal  of  Power  Sources  188  (2009)  170-176 


Table  1 

Reactions  considered  in  this  paper. 


By  applying  above  assumptions  and  simplifications,  the  reduced 
SOFC  model  can  be  developed  based  on  mass  and  energy  balances, 
and  electrochemical  principles.  Fig.  2  gives  the  schematic  represen¬ 
tation  of  the  reduced  model. 


and 

1  dNc  _ 
Ac  dz  ~ 


(5) 


2.1  Mass  balances 


There  are  five  gas  species  in  the  anode  channels,  i.e.,  CH4,  H20, 
H2,  CO  and  C02,  and  two  gas  species  in  the  cathode  channels,  i.e.,  02 
and  N2 .  By  applying  the  molar  balance  to  species  i  in  the  anode  chan¬ 
nels,  the  following  molar  fraction  dynamic  equation  is  obtained: 


d(Ca*u) 

dt 


,  1  3(NaXi,a)_c 
dz 


i  e  (CH4,  H20,  H2,  CO,  C02), 


(1) 


The  molar  source  terms  of  individual  gas  species,  i.e.,  Sm  i  a  and 
Sm,itCf  and  the  total  molar  source  terms,  i.e.,  Sm,a  and  5m)C,  are 
related  to  the  rates  of  reactions  in  the  SOFC.  Defining  a  vector 

j 

r  =  [  q  Hi  Hii  riv  ]  to  represent  the  rates  of  reactions  (i)-(iv), 
these  molar  source  terms  can  be  expressed  as 


W  W  W  W 
A^  Aa  Aa  Aq 


) 


r, 


(6) 


and 


where  Ca  is  the  total  molar  concentration  of  anode  gases,  xi  a  the 
molar  fraction  of  species  i  in  the  anode  gas  stream,  A a  the  total 
section  area  of  anode  channels,  Na  the  molar  flow  rate  of  the  anode 
gas  stream,  and  Sm  i  a  the  molar  source  term  of  species  i  in  anode 
channels. 

For  the  gases  in  the  anode  channels,  the  following  total  molar 
balance  equation  can  also  be  built: 


dCa  1  dNa  = 
dt  +  Aa  dz 


(2) 


where  5m,a  is  the  total  molar  source  term  of  the  gases  in  anode  chan¬ 
nels.  Neglecting  the  total  molar  concentration  transient  of  gases  in 
anode  channels  (i.e.,  the  first  term  on  the  left  hand  side  of  Eq.  (2)), 
we  have  the  following  ordinary  differential  equation: 


1  dNa  = 
A*  dz 


(3) 


Similarly,  for  the  gases  in  the  cathode  channels  we  have: 


d(CcxUc)  1  d{NcxUc) 
dt  Ac  dz 


i  e  {02,  N2}, 


(4) 


_  / w  w  w  w\. 

Vk  =  afcdias(^-^’^-^:Jr’  ks{a,c},  (7) 

where  ai  k  e  R4  is  the  stoichiometric  vector  of  species  i  in  gas  stream 
k  in  reactions  (i)-(iv),  akeR4  the  total  stoichiometric  vector  of  gas 
stream  k  in  reactions  (i)-(iv),  and  W  the  cell  width. 

The  reaction  rates  are  calculated  as  follows  [2,4,10]: 

H  =  krPaXCH4,a  eXP  (^t)  ’  ® 


Hi  —  ^sh^a^CO,a 


^C02,a^H2,aAc0,a^H20^ 


(9) 


Hii  —  riv  —  Tv  —  2 p  ’  (10) 

where  EA  r  is  the  activation  energy  of  the  steam  reforming  reac¬ 
tion,  T  the  SOFC  temperature,  J  the  average  current  density, 
and  Ksh  the  equilibrium  constant  of  the  WGS  reaction  with 
JCsh  =  exp(4276/T-  3.961 )  [11  ]. 


Fuel 


Air 


O  z 


Node 


SOFC 


V7 


Depleted 

fuel 

- ► 

- ► 

Depleted 

air 


Fig.  2.  Schematic  representation  of  the  reduced  model. 
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2.2.  Energy  balance 


According  to  Simplification  (a),  the  energy  balance  equation  for 
the  SOFC  can  be  written  as 

^PENPPENCp,PEN  +  AntPintcp,int  3 T 

A  3 1 

Aa  3(Ca/ia)  Ac  3(Cchc)  1  3 (Na/ia)  1  3 (Nchc) 

X  dt  ~A  3 1  A  3 z  A  3z 


ApEN^PEN  +Aint^int  32T 

A  3z2 


+ 


(11) 


where  A  is  the  cell  section  area,  APEN  (Aint),  pPEN  (pint),  cP)PEN  (cp>int), 
and  Apen  (Aint)  the  section  area,  density,  mass  specific  heat  capacity, 
and  thermal  conductivity  of  the  PEN  (interconnect),  respectively,  ha 
(hc)  the  molar  specific  enthalpy  of  the  anode  (cathode)  gas  mixture, 
and  Sh  the  heat  source  term. 

The  first  term  on  the  left-hand  side  of  Eq.  ( 1 1 )  represents  the  rate 
of  energy  accumulation  in  the  SOLID;  the  second  and  third  terms 
represent  those  in  the  anode  and  cathode  gas  phases,  respectively, 
and  can  be  neglected  because  the  heat  capacity  of  the  gas  phases 
is  far  smaller  than  that  of  the  SOLID;  the  forth  and  fifth  terms  rep¬ 
resent  the  enthalpy  fluxes  due  to  gas  convection  in  the  anode  and 
cathode  channels,  respectively.  The  first  term  on  the  right-hand  side 
of  Eq.  (11)  represents  the  heat  conduction  in  the  SOLID  modeled  by 
Fourier’s  law  of  conduction. 

The  heat  source  term,  Sh,  is  expressed  as 


Sh 


(-n  A -  Hi  AHii  -  rv  A Hv  -UJ)  W 
A 


(12) 


where  A  Hj  is  the  enthalpy  change  of  reaction  j  at  the  temperature  of 
298.15  K,  and  U  the  cell  operating  voltage.  In  Eq.  (12),  Simplification 
(b)  has  been  employed  (i.e.,  the  uniform  current  density  distribu¬ 
tion  is  assumed),  and  thereby  the  electrochemically  generated  heat 
is  also  uniformly  distributed  in  the  SOFC. 

The  specific  enthalpy  of  the  gas  mixture  (i.e.,  fta  and  hc  in  Eq. 
(11 ))  is  calculated  as  a  function  of  the  local  gas  molar  fraction  by 


hk  -  y~Xj  khj,  k e  {a,  c).  (13) 

1=1 


The  specific  enthalpy  of  gas  species  i  at  the  temperature  of  T, 
h/(T),  is  calculated  by 


WT) 


f 

Jjref 


cpj(T )  dT, 


(14) 


with  7^  =  298.15  K,  where  cp  i  is  the  molar  specific  heat  capacity  of 
species  i  and  is  expressed  as  a  third  degree  polynomial  function  of 
T  [12]. 


with  H0(T)  =  1.2723  -  2.7645  x  10_4f,  where  H0(T)  is  the 
reversible  open  circuit  voltage  at  the  temperature  of  T  and 
the  standard  pressure  [11]. 

The  cell  operating  voltage  is  calculated  by 


u  =  Uocv -JRohm(T)  -  —  sinh-1 

/Igf 


2 RT  .  .  _i 

- F  Slnh 

neF 


(16) 


The  second  to  fifth  terms  on  the  right-hand  side  of  Eq.  (16)  rep¬ 
resent  the  ohmic,  anode  and  cathode  activation,  and  concentration 
overpotentials,  respectively.  In  this  equation,  Rohm  is  the  internal 
resistance  of  the  cell,  which  is  a  function  of  the  temperature  and  is 
calculated  from  the  conductivity  of  the  individual  layers  [4] ;  j0,a  and 
jo  c  are  the  anode  and  cathode  exchange  current  densities,  respec¬ 
tively,  and  are  determined  by  the  method  of  Ref.  [13];  jL  is  the 
limiting  current  density  at  which  the  fuel  is  used  up  at  a  rate  equal 
to  its  maximum  supply  speed  in  the  electrode  and  the  cell  voltage 
falls  rapidly  to  zero. 


2.4.  Boundary  conditions 

The  two  ends  of  the  SOLID  are  both  assumed  to  be  insulated.  The 
inlet  boundary  conditions  of  gas  streams  are  specified  by  the  inlet 
gas  conditions  including  the  anode  and  cathode  gas  temperatures, 
molar  fractions  and  flow  rates. 


3.  Solution  method 

Eqs.  (1),  and  (3)-(16),  and  the  relevant  boundary  conditions 
compose  an  equation  system,  which  can  be  solved  using  the  finite 
volume  method  [14].  The  SOFC  domain  is  divided  into  n  uniform 
nodes  along  the  cell  length  direction,  with  the  grid  points  located 
at  the  centers  of  the  nodes.  The  convection  terms  in  Eqs.  (1 )  and  (4) 
are  discretized  using  the  first-order  upwind  scheme,  and  the  dif¬ 
fusion  term  in  Eq.  (11)  using  the  second-order  central  difference 
scheme  [14].  As  the  equation  system  is  stiff,  its  transient  terms 
are  integrated  using  the  stiff  solver  ode23s  provided  by  MATLAB 
R2006a. 

In  each  time  step,  according  to  the  specified  average  current 
density  the  electrochemical  model  is  first  solved  to  obtain  the  cell 
voltage.  Afterwards,  the  molar  source  terms  (Eqs.  (6)  and  (7))  and 
heat  source  term  (Eq.  (12))  are  calculated.  Next,  the  molar  flow 
rates  Na  and  Nc  are  calculated  from  Eqs.  (3)  and  (5).  Then,  the  state 
variables  are  updated  from  Eqs.  (1),  (4)  and  (11).  According  to  the 
solution  strategy  presented,  a  MATLAB  code  is  designed  to  solve 
this  model. 


2.3.  Electrochemical  model 


4.  Model  evaluation 


In  this  model,  the  average  current  density,  which  reflects  the 
external  load  demand,  is  specified,  and  the  cell  operating  volt¬ 
age  is  then  determined.  The  cell  operating  voltage  is  calculated  by 
subtracting  various  irreversible  overpotentials  from  the  reversible 
open  circuit  voltage.  According  to  Simplification  (b),  the  average 
cell  temperature  and  gas  molar  fractions  are  used  in  the  calcula¬ 
tion  of  the  reversible  open  circuit  voltage  and  various  irreversible 
overpotentials. 

The  reversible  open  circuit  voltage  is  given  by  the  Nernst  equa¬ 
tion: 


RT 

Uocv  =  U0(t)+—  In 


X  ^0.5  n0.5 
xH2,a^o2,c1  c 

*H20,a 


(15) 


The  reduced  SOFC  model  (R.M.)  is  evaluated  with  a  detailed 
model  (D.M.)  [15]  and  a  lumped  parameter  model  (L.M.).  The  D.M. 
is  also  a  ID  SOFC  model  as  the  R.M.  but  without  introducing 
Simplifications  (a)  and  (b).  In  the  D.M.,  the  SOFC  is  divided  into 
three  temperature  layers,  i.e.,  the  fuel  stream  in  anode  channels, 
the  air  stream  in  cathode  channels  and  the  SOLID.  The  convec¬ 
tive  heat  transfer  and  the  gas  diffusion  between  the  SOLID  and 
the  gas  channels  are  considered  comprehensively.  Moreover,  the 
Newton-Raphson  iteration  is  employed  to  calculate  the  current 
density  distribution  and  cell  voltage.  The  L.M.  is  a  OD  SOFC  model 
as  given  in  Ref.  [7],  which  neglects  the  spatial  variations  of  SOFC 
variables,  such  as  the  current  density,  temperatures  and  gas  com¬ 
positions,  and  is  built  by  applying  mass  and  energy  balances  to 
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Table  2 

SOFC  specifications. 
Dimensions  of  the  cell 


Cell  length  (m) 

0.1 

Cell  width  (m) 

0.1 

Anode  thickness  (m) 

500  x 10-6 

Cathode  thickness  (m) 

50  x 10-6 

Electrolyte  thickness  (m) 

20  x 10-6 

Interconnect  thickness  (m) 

500  x 10-6 

Fuel  channel  height  (m) 

1.0  x  10-3 

Air  channel  height  (m) 

1.0  x  10-3 

Material  properties  of  the  cell 

Anode  electrical  conductivity  (C2-1  m-1 ) 

8.0  x  104 

Cathode  electrical  conductivity  (£2-1  m-1  ) 

8.0  x  103 

Electrolyte  ionic  conductivity  (£2-1  m-1 ) 

3.34  x  104  exp(- 10,300/7) 

PEN  density  (kg  m-3 ) 

5900 

PEN  heat  capacity  (J  kg-1  K-1 ) 

500 

PEN  thermal  conductivity  (W  irr 1  K-1 ) 

2 

Interconnect  density  (kg  irr3 ) 

8000 

Interconnect  heat  capacity  (J  kg-1  I<-1 ) 

500 

Interconnect  thermal  conductivity  (W  m-1 I<-1 ) 

25 

Activation  and  concentration  overpotential  data 

Activation  energy  of  anode  (J  mol-1 ) 

140,000 

Activation  energy  of  cathode  (J  mol-1 ) 

137,000 

Pre-exponential  coefficient  for  anode  (Am-2) 

0 

0 

X 

0 

06 

Pre-exponential  coefficient  for  cathode  (Am-2) 

1.5  x  1010 

Limiting  current  density  (A  m-2 ) 

12,000 

the  global  SOFC  domain.  The  L.M.  therefore  can  predict  a  sin¬ 
gle  SOFC  temperature.  The  transient  terms  of  the  D.M.  and  the 
L.M.  are  also  integrated  using  ode23s.  A  co-flow  planar  anode- 
supported  DIR-SOFC  is  adopted  in  the  evaluation.  The  specifications 
of  the  SOFC  are  listed  in  Table  2.  All  the  computations  are  per¬ 
formed  on  an  AMD  Sempron™  1.6GFIz  computer  with  512  MB 
RAM. 


The  steady-state  operating  conditions  used  for  model  evalua¬ 
tion  are  listed  in  Table  3,  where  Case  1  is  the  nominal  operating 
condition  and  Cases  2-5  represent  the  operating  conditions  of  high 
load  current,  high  fuel  utilization,  low  air  ratio,  and  high  gas  inlet 
temperature  difference,  respectively.  The  steady-state  performance 
of  the  D.M.  with  100  nodes,  the  L.M.,  and  the  R.M.  with  100  nodes 
under  these  operating  conditions  are  listed  in  Tables  4-6,  respec¬ 
tively.  Among  these  models,  the  D.M.  with  100  nodes  is  considered 
to  be  the  most  accurate  one  because  the  least  simplifications  and 
100  discretized  nodes  are  employed.  The  terms  compared  involve 
anode  and  cathode  outlet  gas  molar  fractions,  molar  flow  rates,  and 
temperatures,  average  and  outlet  SOLID  temperatures,  cell  voltage, 
and  computational  time  per  time  step. 

From  Tables  4  and  5,  the  L.M.  and  the  D.M.  can  first  be  com¬ 
pared  to  show  the  performance  of  the  lumped  model.  It  can  be 
seen  that  the  L.M.  achieves  the  same  prediction  results  as  the  D.M. 
in  the  cathode  outlet  gas  molar  flow  rate  and  molar  fractions.  FIow- 
ever,  there  exist  noticeable  prediction  errors  in  the  anode  outlet  gas 
molar  flow  rate  and  molar  fractions,  particularly  the  outlet  molar 
fractions  of  CFI4,  H2,  and  CO.  This  is  disadvantageous  for  system 
research,  because  in  an  SOFC  system  the  anode  exit  gases  are  usu¬ 
ally  burned  for  energy  recycling.  Moreover,  it  can  be  found  that 
the  temperature  predicted  by  the  L.M.  (e.g.,  1129.3  K  for  Case  1) 
is  close  to  the  outlet  SOLID  temperature  (e.g.,  1121.8  K  for  Case 
1)  rather  than  the  average  SOLID  temperature  (e.g.,  1037.3 1<  for 
Case  1 )  predicted  by  the  D.M.  Therefore,  the  temperature  predicted 
by  the  L.M.  should  be  viewed  as  the  outlet  temperature  for  a  co¬ 
flow  SOFC.  It  is  worth  noting  that  the  prediction  errors  of  Case 
3  in  the  outlet  temperatures  are  all  beyond  30  K  and  that  in  the 
cell  voltage  approximates  0.1  V.  This  indicates  that  the  prediction 
accuracy  of  the  L.M.  significantly  decreases  at  high  fuel  utilization 
level. 


Table  3 

Steady-state  operating  conditions  used  for  model  evaluation. 


Case  number 

1 

2 

3 

4 

5 

Average  current  density  (A  m-2 ) 

5000 

6000 

5000 

5000 

5000 

Fuel  utilization 

0.8 

0.8 

0.9 

0.8 

0.8 

Air  ratio 

7.5 

7.5 

7.5 

5.5 

7.5 

Inlet  fuel  temperature  (K) 

1023 

1023 

1023 

1023 

1048 

Inlet  air  temperature  (I<) 

1023 

1023 

1023 

1023 

998 

Inlet  gas  composition 

Fuel:  CH4  0.33,  H20  0.67;  air:  02  0.21,  N2  0.79 

Operating  pressure  (bar) 

1 

Table  4 

Steady-state  performance  of  the  D.M.  with  100  nodes. 


Case  1 

Case  2 

Case  3 

Case  4 

Case  5 

yOllt 

CH4,a 

6.0400e-5 

1.4616e-4 

6.0926e-6 

4.0926e-5 

1.5841e-4 

YOUt 

H20,a 

0.6772 

0.6785 

0.7398 

0.6793 

0.6767 

vout 

H2,a 

0.1239 

0.1227 

0.0614 

0.1219 

0.1244 

vOUt 

XCO,a 

0.0349 

0.0358 

0.0181 

0.0370 

0.0340 

vOUt 

C02,a 

0.1639 

0.1629 

0.1807 

0.1618 

0.1647 

N°ut  (mols-1) 

4.0726e-4 

4.8863e-4 

3.6205e-4 

4.0728e-4 

4.0718e-4 

vout 

02,c 

0.1872 

0.1872 

0.1872 

0.1786 

0.1872 

vOUt 

N2,c 

0.8128 

0.8128 

0.8128 

0.8214 

0.8128 

N°ut  (mols-1) 

4.4974e-3 

5.3968e-3 

4.4974e-3 

3.2635e-3 

4.4974e-3 

fs  (I<) 

1037.3 

1043.8 

1053.3 

1037.4 

1021.1 

Tsout  (K) 

1121.8 

1135.2 

1133.9 

1149.2 

1112.0 

Tout  (K) 

1122.5 

1136.0 

1134.4 

1150.0 

1112.7 

Tcout  (K) 

1116.6 

1127.0 

1130.3 

1145.0 

1106.0 

U(V) 

0.6582 

0.6238 

0.6599 

0.6626 

0.6252 

Ar  (s) 

23.6102 
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Table  5 

Steady-state  performance  of  the  L.M. 


Case  1 

Case  2 

Case  3 

Case  4 

Case  5 

o  u 

0.0111 

0.0123 

0.0077 

0.0094 

0.0123 

y  0  Lit 

H20,a 

0.7046 

0.7081 

0.7603 

0.7019 

0.7067 

X°Llt 

H2  ,a 

0.0922 

0.0882 

0.0379 

0.0956 

0.0896 

vout 

*CO,a 

0.0259 

0.0254 

0.0121 

0.0287 

0.0241 

vOUt 

XC02,a 

0.1662 

0.1660 

0.1820 

0.1644 

0.1673 

N°ut  (mols-1) 

3.9845e-4 

4.7700e-4 

3.5657e-4 

3.9977e-4 

3.9752e-4 

vout 

02,c 

0.1872 

0.1872 

0.1872 

0.1786 

0.1872 

yOUt 

*N2,c 

0.8128 

0.8128 

0.8128 

0.8214 

0.8128 

N°ut  (mols-1) 

4.4974e-3 

5.3968e-3 

4.4974e-3 

3.2635e-3 

4.4974e-3 

Ts  (I<) 

1129.3 

1138.7 

1165.1 

1152.1 

1115.5 

Tsout  (K) 

1129.3  (7.5)a 

1138.7(3.5) 

1165.1  (31.2) 

1152.1  (2.9) 

1115.5  (3.5) 

Taout  (K) 

1129.3  (6.8) 

1138.7(2.7) 

1165.1  (30.7) 

1152.1  (2.1) 

1115.5(2.8) 

Tcout  (I<) 

1129.3(12.7) 

1138.7(11.7) 

1165.1  (34.8) 

1152.1  (7.1) 

1115.5(9.5) 

U(V) 

0.6405  (-0.0177) 

0.6124  (-0.0114) 

0.5632  (-0.096) 

0.6648  (0.0022) 

0.6200  (-0.0052) 

Ar  (s) 

0.0137 

a  The  value  in  parenthesis  is  the  prediction  error  compared  with  the  D.M. 


Table  6 

Steady-state  performance  of  the  R.M.  with  100  nodes. 


Case  1 

Case  2 

Case  3 

Case  4 

Case  5 

yOllt 

ACH4,a 

3.7721e-5 

9.3690e-5 

6.6218e-6 

1.6553e-5 

8.7239e-5 

y  °  Lit 

H20,a 

0.6767 

0.6777 

0.7396 

0.6785 

0.6760 

vOUt 

XH2,a 

0.1245 

0.1235 

0.0616 

0.1227 

0.1252 

yOUt 

“TO,  a 

0.0344 

0.0352 

0.0179 

0.0363 

0.0336 

vOUt 

XC02,a 

0.1644 

0.1635 

0.1809 

0.1625 

0.1652 

N°ut  (mols-1) 

4.0728e-4 

4.8868e-4 

3.6205e-4 

4.0730e-4 

4.0724e-4 

vOUt 

X02,c 

0.1872 

0.1872 

0.1872 

0.1786 

0.1872 

yOUt 

XN2,c 

0.8128 

0.8128 

0.8128 

0.8214 

0.8128 

N°ut  (mols-1) 

4.4974e-3 

5.3968e-3 

4.4974e-3 

3.2635e-3 

4.4974e-3 

Ts  (I<) 

1041.1  (3.8) 

1047.4  (3.6) 

1050.0  (-3.3) 

1049.0(11.6) 

1028.1  (7.0) 

Tsout  (K) 

1116.3  (-5.5) 

1127.0  (-8.2) 

1128.2  (-5.7) 

1141.0  (-8.2) 

1105.6  (-6.4) 

Taout  (I<) 

1116.3  (-6.2) 

1127.0  (-9.0) 

1128.2  (-6.2) 

1141.0  (-9.0) 

1105.6  (-7.1) 

Tcout  (I<) 

1116.3  (-0.3) 

1127.0(0.0) 

1128.2  (-2.1) 

1141.0  (-4.0) 

1105.6  (-0.4) 

U(V) 

0.6614(0.0032) 

0.6268  (0.0030) 

0.6678  (0.0079) 

0.6740  (0.0114) 

0.6287  (0.0035) 

Ar  (s) 

3.4396 

Then,  from  Tables  4  and  6,  the  R.M.  and  the  D.M.  can  be  com¬ 
pared  to  show  the  performance  of  the  reduced  model.  It  can  be 
found  that  the  R.M.  also  achieves  the  same  prediction  results  as 
the  D.M.  in  the  cathode  outlet  gas  molar  flow  rate  and  molar  frac¬ 
tions.  Moreover,  there  is  very  little  difference  in  the  anode  outlet 
gas  molar  flow  rate  and  molar  fractions  predicted  by  the  two  mod¬ 
els.  This  is  an  improvement  compared  with  the  L.M.  The  outlet 
temperatures  predicted  by  the  R.M.  are  generally  lower  than  those 
predicted  by  the  D.M.  This  is  because  introducing  the  simplifica¬ 
tion  of  uniform  current  density  distribution  flats  the  temperature 
distributions.  As  a  whole,  the  prediction  errors  in  the  average  and 
outlet  temperatures  and  the  cell  voltage  under  different  operat¬ 
ing  conditions  are  modest.  The  prediction  errors  in  these  variables 
slightly  increase  in  the  low  air  ratio  case  (i.e.,  Case  4)  and  the  high 
gas  inlet  temperature  difference  case  (i.e.,  Case  5).  These  results 
indicate  that  the  R.M.’s  accuracy  is  not  significantly  compromised 
by  the  two  simplifications  under  the  operating  conditions  investi¬ 
gated.  It  can  also  be  found  that  the  computational  time  of  the  R.M.  is 
significantly  lower  than  that  of  the  D.M.  due  to  employing  the  two 
simplifications. 

To  further  verify  the  dynamic  performance  of  the  reduced 
model,  Fig.  3  shows  the  dynamic  responses  of  the  D.M.  with 
100  nodes,  the  L.M.,  and  the  R.M.  with  100  nodes  to  the 
applied  average  current  density  disturbances.  It  can  be  found 
that  these  models  exhibit  similar  dynamics  in  the  outlet  SOLID 


Fig.  3.  Comparison  of  dynamic  responses  of  the  D.M.  with  100  nodes,  the  L.M.,  and 
the  R.M.  with  100  nodes  to  applied  load  current  disturbances. 
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temperature  and  the  cell  voltage  for  the  applied  load  current  dis¬ 
turbances. 

When  the  node  number  of  the  R.M.  is  reduced  to  one,  it  will 
have  the  same  equation  form  as  the  L.M.  Therefore,  the  L.M. 
can  be  viewed  as  the  R.M.  with  one  node.  The  R.M.  approxi¬ 
mates  the  spatially  distributed  nature  of  SOFCs  along  the  gas 
flow  direction  by  connecting  a  series  of  scale-down  L.M.s  with 
consideration  of  the  heat  conduction  and  gas  transport  between 
adjacent  nodes.  Therefore,  the  node  number  can  be  taken  as  an 
adjustable  parameter  in  modeling,  and  is  determined  according 
to  specific  modeling  requirements  for  accuracy  and  computational 
time. 

5.  Conclusions 

A  reduced  ID  dynamic  model  of  a  planar  DIR-SOFC  has  been 
presented  in  this  paper  for  system  research  by  introducing  two  sim¬ 
plifications.  A  detailed  1 D  SOFC  model  without  introducing  the  two 
simplifications  and  a  lumped  parameter  SOFC  model  are  employed 
to  evaluate  the  present  model  for  a  co-flow  planar  anode-supported 
DIR-SOFC.  From  the  model  evaluation  results,  several  conclusions 
can  be  drawn  as  follows: 

(1)  Under  the  operating  conditions  investigated,  the  accuracy  of 
the  reduced  model  is  not  significantly  compromised  by  the 
two  simplifications  in  prediction  of  the  outlet  gas  flow  rates 
and  molar  fractions,  the  outlet  temperatures,  and  the  cell  volt¬ 
age,  while  its  computational  time  is  significantly  decreased  by 
them.  More  noticeable  predictions  errors  are  possible  under 
extreme  operating  conditions.  Moreover,  it  is  quite  simple  in 
form.  Therefore,  the  reduced  model  provides  a  choice  of  mod¬ 
eling  for  SOFC  system  research. 

(2)  Compared  with  the  lumped  model,  the  reduced  model  is  an 
improvement  with  regard  to  accuracy,  because  the  reduced 
model  can  be  viewed  as  a  series  of  scale-down  lumped  models, 


which  are  connected  to  approximate  the  spatially  distributed 
nature  of  SOFCs  along  the  gas  flow  direction. 

(3)  The  discretized  node  number  for  solving  the  reduced  model 
can  be  taken  as  an  adjustable  parameter  in  modeling,  and  is 
determined  according  to  specific  modeling  requirements  for 
accuracy  and  computational  time. 

In  the  future,  the  reduced  SOFC  model  will  be  used  for  system 
optimization  and  control  design. 
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